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Abstract 

The spin and charge structure factors are calculated for the Hubbard model on 

the square lattice near half-filling using a spin-rotation invariant six-slave boson 

representation. The charge structure factor shows a broad maximum at the zone 

corner and is found to decrease monotonically with increasing interaction strength 

and electron density and increasing temperature. The spin structure factor develops 

with increasing interaction two incommensurate peaks at the zone boundary and 

along the zone diagonal. Comparison with results of Quantum Monte Carlo and 

variational calculations is carried out and the agreement is found to be good. The 

limitations of an RPA-type approach are pointed out. 
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I. INTRODUCTION 



Soon after the discovery of superconductivity in the cuprate materials, it was suggested 
[ ij that this phenomenon is closely related to strong correlation effects. Indeed correlations 
are responsible for the insulating state observed in the parent compounds. The simplest 
Hamiltonian accounting for such Mott insulators is the one band Hubbard model. It poses a 
serious challenge to the theoretician since ordinary many-body perturbation theory breaks 
down for strong coupling, being unable to account for Mott insulator state. A number of 
new techniques have been developed, either fully numerical such as Quantum Monte Carlo 
calculations or exact diagonalizations of small systems [§, or analytical using the Hubbard 
X-operator technique (for a recent work see |§), the self-consistent 2-particle theory [|J, 
the dynamical mean field approximation || or slave bosons. The slave boson method has 
been applied to a whole range of problems with local Coulomb interaction: the Kondo 
impurity model H0], the Kondo lattice model |7|-[10|, the Anderson Hamiltonian |6|,|TT| the 
Hubbard model fll2| , [T3H possibly with orbital degeneracy |14]] and even the Bose-Hubbard 
model fll5f . In the Kotliar and Ruckenstein (KR) slave boson technique |12 the Gutzwiller 



Approximation p6Hl9|] appears as a saddle-point approximation of this field theoretical 
representation of the Hubbard model. In the latter a metal-insulator transition occurs at 
half-filling in the paramagnetic phase as discussed by Lavagna PC[ ]. The contribution of 
the thermal fluctuations has been calculated f2]| and turned out to be incomplete as this 
representation, even though exact, is not manifestly spin-rotation invariant. Spin-rotation 
invariant [22] and spin and charge-rotation invariant |23f formulations have been proposed, 
all sharing the advantage of treating all the atomic states on an equal footing, and the 



first one was used to calculate correlation functions [24| and spin fluctuation contributions 



to the specific heat p5| . Comparisons of ground state energy with Quantum Monte-Carlo 
simulations, including antiferromagnetic ordering |26| and spiral states p7f , or with exact 
diagonalization data have been done and yield excellent agreement, and a magnetic phase 
diagram has been proposed |^9]|!0|. The magnetic susceptibility has been evaluated ]3l|] and 
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shows in the strong coupling regime a maximum in its doping dependence in agreement 
with the observed behavior in the cuprates. Such a field theoretical description is especially 
useful since it allows for the calculation of dynamical quantities as well. This is the goal of 
the paper. We first derive expressions for the spin and charge auto-correlation functions, 
which we then evaluate numerically and determine the structure factors. We discuss how 
they depend on the doping, the interaction strength and the temperature. 



II. RESPONSE FUNCTIONS 

In this work we calculate the spin and the charge structure factors of the Hubbard model 
on the square lattice within the Spin Rotation Invariant (SRI) slave boson formulation 
23| . After having compared ground state energies, effective band widths and dispersions 



and magnetic susceptibilities with results obtained by other techniques, generally with very 
good agreement, structure factors provide another test for the approach, which is free of any 
adjustable parameter. 

In the Spin Rotation Invariant slave boson representation, the Hubbard Hamiltonian is 
expressed in terms of slave boson operators e^, p. dj for empty, singly occupied and doubly 



occupied sites and pseudo-fermions operators f a as |12|j2^ , [23| 

# = E tu £ *i»j£f s j*jw + U H dfd> (i) 

where 

Zi = efLM^ P t + PtLMEidi (2) 

with 

Mi = [1 + (efei + pf Q pio + pt ■ Pi + dfdi)}^ 

Li = [(l-dM)ro-P^r 1/2 

Ri = [(l-e^)r -^l-V2 (3) 



3 



The under-bar denotes a 2 x 2 matrix in spin space, r is the unit matrix and p. is the time 
reverse of operator p.. As usual the slave boson operators have to fulfill constraints. Here 
they read 

e t e i + (pfoPio +Pi -Pi) + dfdi = 1 

fl (J 

(pfoPi+PiPiO-iPi ><Pi) = Y.^a'fZ'f™ ( 4 ) 



la 

i 

crcr 



and they are respectively enforced by the constraint fields a,/3 ,(3. The scalar and vector 
slave boson fields pio and pi are defined as 

Pi = ^(piolv + PiZ) (5) 

where r is the vector of Pauli matrices. The paramagnetic mean-field free energy results 
into 

F = -2T J dep(e) ln[l + exp —E/T) + Ud 2 

+a(e 2 + d 2 +p 2 -l)- (3 (p + 2d 2 ) . (6) 

After having fixed the notation we can proceed to the calculation of the dynamical suscep- 
tibilities. The linear response to an external field is given by the one-loop order calculation. 
The expression for the dynamical spin susceptibility has already been obtained by Li et al 



24j| , but the calculation of charge fluctuations turned out to be more involved. In early 
calculations it has been assumed that the gauge symmetry group of both KR and SRI 
representations allows for gauging away the phases of all slave boson fields |2l| , |20| . It later 
turned out that this conclusion is erroneous, and one Bose field has to be complex p2| , |l3| , p3| . 
Choosing it to be the field describing double occupancy allows for describing the physics of 



the upper Hubbard band p3| . We note that the fully symmetric gauge approach has been 
applied to the calculation of the Landau parameters for liquid 3 He by Li and Benard p4 |. 
The spin and charge auto-correlation functions Xs and \ c can be obtained out of: 

^(i-iV-r')=<f[n iff (T)n /y (r')] > (7) 



as 

Xs(i-i',r -t) = ^oo B cc ,{i-i ,t -t) 

t 

(7(7 

Xc(i-i',T-T) = ^^(i-^T-r') (8) 

era 

Using the constraints and the mapping 

1 3 1 

n i° = 2 J2ptPi» + d t d i + ^PmP^ + PtsPio ~ i(PaPa ~ PaPn)] , (9) 

li=0 



one can express the density fluctuations in terms of the slave boson fields as 
]T5n CT = 5{d+d-e+e) = 5N 

a 

J2 aSn <r = $(PoP3 +P3P0) =5S . (10) 

The correlation functions can be written in terms of the slave boson correlation functions 
as: 

Xs(k) = E CT y oa < 5n a (-k)5n a ,(k) > =< 5S(-k)5S(k) > 

Xc(k) = Y. aa > < 6n a (-k)5n a , (k) > =< 5N(-k)5N(k) > . (11) 

Performing the calculation to one-loop order, one can make use of the propagators given in 
the appendix to obtain: 

X c(k) = 2e 2 S n \k) - AedS u \k) + 2d 2 S^(k) 

Xs (k) = 2p 2 S^(k) . (12) 
Including the inverse matrix elements, we get 



and 



v ^3 3 (A;) 1 S 55 (A;)e 2 [-2^r 1 (A;) + 8d Po T 2 (k) - 8d 2 T 3 (k)} - 2e A p 2 S 2 5 (k)u 2 ^ 

uf, ' ) : ~ Snikxr^r^k) - ri(fc)) + ^r 3 (k)s 55 (ky " (L3) 



XA ) 1 + A iXo (k) + A lX i(k) + A 2 [ x i(k) - Xo(k) X 2(k)} ' 1 ] 
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where 



At 
A 2 



(Zpg)" 1 



a - Pq + e z 



«9 2 z T 
dpi 



.<9z 



TN2 



-i 

P Z OT— 

ops 



dz 



T\2 



(15) 



and 



iMfc) = -S 55 (A;)(e 2 S 22 (A ; )-2e^ 12 (A;) + C / 2 5 11 (A;)) + (eS 2 5(A ; )-^ 1 5(A;)) 2 

r 2 (fe) = -5 , 55 (A;)(e 2 S , 2 4(A;) -poeSi2(fc) - edSu(k) + dp S n (k)) 

+ (eS , 25 (/;;) - rfS , i 5 (A;))(eS , 45(A;) -pd-SisW) 

T 3 (A;) = -5 55 (fc)(e 2 5 4 4(A;)- 2e Po S u (k) + plS n (k)) + (eS 45 {k) - p S 15 {k)) 2 . (16) 

In the following we shall evaluate numerically the spin and charge structure factors: 
dw ImXxiQi u + £0) 



a; = s, c 



(17) 



7T 1 — exp(—uj/T) ' 

Note that while deriving the expressions (|T2| ) for the spin and charge correlation func- 
tions, we are dealing with both complex (d) and real (e,p^) fields. For the latter, following 
an argument by Read and Newns 0, the contribution from the measure to the action 
{J2i, n m ( e i,n YluPin)) nas been neglected. This aspect has been recently re-investigated in 
the framework of the 1/N expansion of the large U Hubbard model p5| . It has been shown 
that including this contribution only leads to a minor change of the action, leaving its 
numerical value unchanged. This argument holds in the present context as well. 



III. RESULTS 

We first determine the spin structure factor as given by Eq. fllTI). In the weak coupling 
limit our result reduces to the RPA |[36|| , and thus for small coupling and particle densities, 
the agreement with the exact solution is expected to be very good. However the RPA is 
getting less and less reliable as the interaction strength increases. This leads to unreasonable 
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results in RPA as e.g. a magnetic instability of the paramagnetic state for any particle filling 



above a critical coupling [ T7| . Such a deficiency is corrected in slave boson mean field theory. 



At zero temperature magnetic instabilities appear only beyond a certain particle filling J3D[ . 
For increasing temperature the region of magnetic long-range order is shrinking rapidly 
pifl . Strictly speaking, long-range order is absent at any finite temperature due to thermal 
excitations of spin waves (Mermin- Wagner theorem), which are not taken into account in 
slave boson mean-field theory. 

Here we display results in the paramagnetic state, away from the instability line. Let us 
start by comparing (Fig. 1) our result for the spin structure factor to Quantum Monte Carlo 
results ||38|| . Here the structure factor is plotted as a function of wave number along straight 



lines T-X-M-T in the Brillouin zone. For U = At, temperature T = t/6 (J3 = 6; (3 — t/T) 
and a doping 8 = 0.275, the overall agreement is good and the trends of the simulations 
are reproduced. Especially the peak is very broad and is centered around (37r/4,7r) for 
the QMC data. Out of our approach we get two peaks located away from (jr, 7r), which 
are indicating the onset of incommensurate short ranged spin order. As compared to the 
QMC calculations we obtain a second peak along the diagonal of the Brillouin zone which 
cannot be resolved in the simulations due to the small size of the system. The height 
of the peak which is located along the zone boundary is larger than the one of the peak 
which is lying on the diagonal of the Brillouin zone. We note that we are comparing the 
slave boson results for an infinite system to the raw QMC data on the 8x8 lattice. We 
believe that this comparison is meaningful since we observed that finite size effects are 
small (in the percent range as compared to the 6x6 lattice) for the set of parameters we 
are using. The agreement between the 2 approaches is mostly qualitative, and the slave 
boson calculation tends to overestimate the tendency towards magnetic ordering. But it 
also allows for gaining additional informations which are not revealed by the simulations 
like the presence of a second peak along the diagonal of the Brillouin zone. QMC data are 
not available in the complete parameter range. As an alternative there exist variational 
approximations at T = to which we can compare our data. A very promising variational 



ansatz is provided by a generalized Baeriswil-Gutzwiller wave- function [39]. In Fig. 2 we 



compare our result for U — It and U — At {(3 — 8, because of onset of long range order in 



our mean-field treatment) and 5 = 0.218 with the zero temperature results of ref. | 3"9]| , and 
we reach the same conclusions as when comparing with the QMC simulations. On top of 
that we notice that the agreement is better for weaker interaction. Note that the position of 
the peaks is temperature dependent. At low temperature it is systematically located away 
from (it, it), except at half-filling, but it moves towards (it, 7r) for increasing temperature, 
as shown on Fig. 3 for U = At. The first result of a small increase of the temperature is 
to suppress the value of the structure factor at its peak position, while the one at the zone 
corner is increasing, up to the point where it becomes the dominant one, as emphasized 
in the inset of Fig. 3. Increasing the temperature further results into an overall reduction 
of the structure factor for large momenta, and an overall increase of it for small momenta. 
The fluctuation-dissipation theorem relates the q = value of the spin structure factor to 
the magnetic susceptibility. As shown in the inset of Fig. 3, the latter is growing with 
increasing temperature for low T, reaches a maximum at (3 ~ 4, and decreases beyond. This 
behavior is reminiscent of the spin gap behavior observed in the High T c superconductors, 
although there it occurs at a smaller energy scale. As one can see from Fig. 3 , the spin 
structure factor appears to approach its zero temperature limiting value at (3 ~ 8. We 
may therefore use the result for U = At at (3 = 8 instead of the T = (j3 = oo) result 
(which is not accessible because of onset of long range order in our mean-field treatment), to 
compare with the variational Monte Carlo (VMC) result at doping concentration 5 = 0.218 
(Fig. 2). Whether incommensurate long range order occurs off half-filling using Otsuka's 



wave-function is not known. In any case it sets in using the Gutzwiller wave-function |40[ . 

We now consider the dependence of the spin structure factor on the interaction and the 
density. Here we fix the temperature to (3 — 8, and we calculate S a (q) for U = t and 
U = 2t and display it on Fig. 4 for S = 0.1, on Fig. 5 for 5 = 0.2 and U/t = 1, 2, 3 
and 4, and on Fig. 6 for 5 = 0.3 and U/t = 2, 4, 6, 7 and 8. In all cases we obtain that 
raising up the interaction generates more spin ordering, and the structures in the curves are 
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much more pronounced. At this temperature the position of the peak remains at the zone 
corner for small doping (see Fig. 4), while it is shifted away from its commensurate value 
for larger doping (see Fig. 5 and 6). The influence of the interaction is stronger when the 
system is denser, and this effect is enhanced by the vicinity of the perfect nesting point. 
For large doping the peaks become much broader resulting in an incommensurate very short 
ranged spin order. Also increasing the hole doping shifts the positions of the maxima of 
S s (q) further away from (71", 7r), both along the M — X line and along the V — M line. One 
may ask to what extent the results presented above can be obtained in the framework of 
an RPA-type (or Fermi liquid type) scheme. For this purpose we compare the result for Xs 
obtained by replacing the effective interaction AAq) by its q = limit and putting A\ and 



A 2 equal to zero in Eq. ( ]T4]JT5|) . This allows for investigating the influence of the momentum 
dependence of the effective interaction on the spin order. The result is displayed on Fig. 
7 for U/t = 2,4 and 6. Clearly the approximation is very good for moderate couplings, 
but gets gradually worse for increasing interaction. Thus the dispersion of the effective 
interaction has a negligible influence on the structure factor for moderate interaction, but 
an important one for intermediate to large couplings, where it strongly shifts and suppresses 
the peaks of the structure factor. We thus conclude that knowing the Landau parameter 
Fq is not all what is needed in order to determine the spin structure factor, especially when 
the fluctuations have a dominant short wavelength character. 

We now turn to the charge structure factor. In Fig. 8 we compare our result for U = 4t 
and S = 0.275 at temperature T = t/6 with the Quantum Monte Carlo result of Dzierzawa 
f38fl . The charge structure factor consists of one broad peak which is centered at (vr,7r). 
Fig. 8 clearly shows that the agreement between both approaches is excellent, and that 
the difference does not exceed a few percent. In Fig. 9 our result at zero temperature and 
5 = 0.218 for U = At, U = 8t, and U = 16t are compared with the variational Monte 



Carlo result of |39[, and the dependence on U is also displayed on Fig. 10. As expected, 
increasing the interaction strength U leads to a suppression and a further broadening of the 
charge structure factor. In Fig. 11, the dependence on doping is shown for U — 8t at zero 
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temperature. In weak coupling, one would expect \ c to decrease upon doping. However the 
opposite behavior holds in a dense strongly correlated system, and our approach succeeds 
in obtaining this subtle effect. As compared to the 2-particle self-consistent theory Q], the 
reason of the success is however quite different. Indeed if we use the "RPA" approximation 
(the same procedure as above, but now in the charge fluctuation sector), we again obtain that 
the difference as compared to the full expression grows under an increase in the interaction 
strength as displayed in Fig. 12. Thus in our theory it is essential to take the dispersion 
of the interaction into account in order to obtain a good result, while such a dispersion is 
neglected in the 2-particle self-consistent theory ||. We also performed the calculation at 
finite temperature. In contrast to the spin structure factor, the charge structure factor is 
mostly temperature independent at low T. This is a pure interaction effect, since in the non- 
interacting limit both structure factors have the same T-dependence. A further increase in 
T simply leads to an overall reduction of the charge structure factor for large momenta, and 
an increase of it for small momenta, as shown in Fig. 13, where S c (q) is plotted for U = At 
and S = 0.275 at temperatures ranging from T = t/1.3 down to T = t/8. Furthermore the 
temperature has no influence on the position of the peak of the charge structure factor. The 
fluctuation-dissipation theorem relates the charge structure factor at q = to the charge 
susceptibility. As indicated in the inset of Fig. 13, the latter decreases monotonously with 
increasing Temperature, in agreement with Fermi liquid theory. 



IV. SUMMARY 

In this work we derived and evaluated the spin and charge structure factors of the Hub- 
bard Model within the Spin Rotation Invariant six-slave boson formulation of the Hubbard 
model. We considered Gaussian fluctuations about the paramagnetic saddle-point, at small 
hole doping and for finite temperature, where the mean field solution is the paramagnetic 
one. The agreement with available exact numerical results for finite size systems and for 
variational wave functions was found to be very good. It is found that increasing either 
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the interaction strength or the density leads to more incommensurate spin order, and less 
commensurate charge order. The slave boson mean field theory appears to provide a good 
starting point for describing dynamical correlations on the level of Gaussian fluctuations. 
It would be of interest to extend the calculations to the magnetically ordered phases, dis- 
cussed e.g. in ref. [28]. Also the consideration of short range antiferromagnetic order [4~1] is 



presumably important at low temperature, and should be taken into account. 
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VI. APPENDIX 



In second order in the bosonic variables the action is given by ||23|| : 

s= E M-q)^Ai)M<i) (is) 



with = e,ip 2 = d',ip 3 = d",ijj 4 = p ,ip 5 = f3 ,ip G = a,ip 7 = pi,^s = Pi, ^9 = P2,^io = 
/3 2 ,ipn = J>3?V ; i2 — /3 3 - Here d! and d" are the real and imaginary parts of the complex 
d-field. The propagator matrix decouples into 4 blocks, one for the charge fluctuations and 
3 for the spin fluctuations. Since the calculation is straightforward and most results can be 
gathered from ref ||24||, we only quote the results. The charge part of the fluctuation matrix 



is given by: 



S n (k) = a + Sn(k) 

S 22 {k) = a - 2f3 + U + S 22 (k) 

S 33 (k) = a-2(3 + U + S' 33 (k) 
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S 23 (k) 




Su(k) 


= a - (3q + S M (k) 


Sij(k) 


= Sij(ky,i^j;i,j = 1,2,4 


S 15 (k) 


1 dz 


S 16 (k) 


= e 


S 25 (k) 


= - 2d -2 Xl{k) dd/ Z ° 


S 26 (k) 


= d 


S 45 (k) 


1 , , . dz 


S* 6 (k) 


= Po 


S 5 5(k) 


= ~\xo(k) 



and the spin part by: 

S 77 (k) = S 99 (k) = Su.u(k) = a-[3 + S 77 (k) 



S 88 (k) = S w>w (k) = S 12>12 (k) = ~Xo(k) 



(19) 



S 7 s(k) = S 9:10 (k) = S u ,i2(k) = ~Po ~ ^i(^)^o • (20) 

We also defined: 

~ 0^ z ! dz dz 

/ , \ cP 1 z r 1 o ' / 1 \ 1 0z^~ dz , „ „ , 

s » w = F -° z °WaJ> + le * ~ 2 z ° x ><- k)] a¥W (21) 

and 

H = EV-^o(p) (22) 

per 

in terms of the pseudo-fermion Greens function Go(p) = \/{iu n — E%) and the effective 
dispersion: 

E n =$t n -n + Po ■ (23) 
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We also introduced the dynamical response functions of the fermionic system: 



Xn(k) = -T,( t p + Wk) nG oAp)G Ap + k) = X n(-kKn = 0,1,2) 

per 

X 2 {k) = -E(W-V) 2 G 0ct (p)G 0ct (p + A;) • (24) 



per 



In Eq. the derivatives are given by: 

dz /- , 2xd , 
V2 Po r](l + 



7/ 



iV2p f] 



dd 

d 2 z 2V2p oV , n . ^3xd\ 

ddP = -TTT {2d + x + —5< 

d 2 z ,2y/2p dr] 



dd'dd" 1 + 5 

d 2 z 2\/2p xr] 

'dF 2 ~ 1 + 5 
d 2 z „ r- , e d 



2V2p r](- + + 2edxrj 2 ) 



dd'de ^ u,v l-5 1 + 5 

()1 ' ~ V2 V (1 + 2 P y + ^- + + 2p 2 xdrf) 



dd'dp ,v " u ' 1 + 5 (1 + 5)' 

d 2 z .2\Z2p Q er] 
dd"de = 1 1-5 
d 2 z 

iV2r)(l + 2p 2 r] 2 ) (25) 



dd"dp 

With x = e + d and rf = 1/(1 — 5 2 ). The other derivatives can be found in Ref. pi ]. This 
extends the result of Bang et a/ [0] to the metallic regime. 
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FIGURE CAPTIONS 



FIG. I. Comparison of the Quantum Monte Carlo (circles) and Slave Boson (full line) 
spin structure factors for U = At, 5 = 0.275 and (3 — 6. 

FIG. 2. Comparison of the Variational Monte Carlo (circles for U — 2t and squares for 
U = At) and Slave Boson (dashed line for U = 2t, (3 = 8 and full line for U = At, (3 = 8) 
spin structure factors for 5 = 0.218 and T = 0. 

FIG. 3. Temperature dependence of the spin structure factor for U = At, and 5 = 0.275. 
The temperatures are (3 = 8 (full line), 4 (dotted line), 2 (dashed-dotted line) and 1.3 (long 
dashed- short dashed line) . Left inset: Magnification of the latter around the M-point. Right 
inset: Temperature dependence of the magnetic susceptibility for U — At and 5 = 0.275. 

FIG. 4. U-dependence of the Slave Boson spin structure factor for U = It (circles) and 
U = 2t (triangles), for 5 = 0.1 and (3 = 8. 

FIG. 5. U-dependence of the Slave Boson spin structure factor for U = It (circles), 
U = 2t (triangles), U = 3t (pluses) and U = At (x's), for 5 = 0.2 and (3 = 8. 

FIG. 6. U-dependence of the Slave Boson spin structure factor for U — 2t (circles), 
U — At (triangles), U = 6t (pluses), U — 7t (x's) and U = 8t (diamonds), for 5 = 0.3 and 
(3 = 8. 

FIG. 7. Comparison of the Slave Boson (x's resp. diamonds resp. stars) and RPA (circles 
resp. triangles resp. pluses) spin structure factors for U — 2t resp. U = At resp. U = 6t, 
and 5 = 0.3 and (3 = 8. 

FIG. 8. Comparison of the Quantum Monte Carlo (triangles) and Slave Boson (full line) 
charge structure factors for U = At, 5 = 0.275 and (3 = 6. 

FIG. 9. Comparison of the Variational Monte Carlo (circles, pluses, diamonds) and 
Slave Boson(full, dashed, dashed-dotted line) charge structure factors for U = At, U = 8t 
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and U = 16*, and 5 = 0.218 and T = 0. 

FIG. 10. U-dependence of the Slave Boson charge structure factor for U = 2t (circles), 
U = 4t (triangles), U = 8t (pluses) and U = 16t (x's), for 5 = 0.275 and T = 0. 

FIG. 11. Doping-dependence of the Slave Boson charge structure factor for 5 = 0.1 
(circles), 5 = 0.2 (triangles), and 5 = 0.3 (pluses), and for U = St and T = 0. 

FIG. 12. Comparison of the Slave Boson (circles, pluses, diamonds) and RPA (squares, 
triangles, x's) charge structure factors for U — 2t resp. U = At, resp. U = 8t, and 5 = 0.218 
and T = 0. 

FIG. 13. Temperature dependence of the charge structure factor for U = At, and 5 = 
0.275. The temperatures are (3 = 8 (circles), 6 (triangles), 4 (pluses), 2 (x's) and 1.3 
(diamonds). Inset: Temperature dependence of the charge susceptibility for U = At and 
5 = 0.275. 
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